1 Introduction

Previously, we observed that time at room temperature (RT) prior to cryopreservation is a source of technical artifacts in blood-based single-cell RNA-seq (scRNA-seq) data. In this suproject, we aim to assess if culturing PBMCs after thawing them and before sequencing can remove the cold-shock signature.

2 Description of the data

To that end, we drew blood from four healthy donors (3 males and 1 female) and cryopreserved the samples after 0h (fresh), 8h or 24h at RT. Subsequently, we thawed the samples and either processed them immediately with 10X Chromium (day 0), or after 2 days of cell culture (day 2). To eliminate batch effects, detect doublets and reduce the library cost, we performed cell hashing. In this protocol, each condition (in our case time-points) is labeled with a specific hashtag oligonucleotide (HTO) that is crosslinked with an antibody. The antibodies bind to ubiquitous cell surface markers, and the HTO are sequenced alongside the single-cell gene expression libraries. We have the following libraries:

  • Day 0: multiplexed with the two donors (male/female1 OR female2/female3), and the 3 time-points per donor (0, 8, 24h).
  • Day 1/2: same as before but after 1 (BCLLATLAS_14 experiment) or 2 (BCLLATLAS_11) days in culture. Moreover, T cells were activated with anti-CD3 antibodies.

2.0.1 Objective

The objective of this notebook is to demultiplex the barcodes (cells) back to its original time-point. To achive that, we will follow the pipeline from Seurat.

3 Pre-processing

3.1 Load packages

library(scater)
library(SingleCellExperiment)
library(Seurat)
library(Matrix)
library(tidyverse)

4 Demultiplex

# Load expression matrix, gene and cell metadata
libraries <- c(
  "BCLLATLAS_11/Tcell_activation_day0", 
  "BCLLATLAS_11/Tcell_activation_day2", 
  "BCLLATLAS_14/Tcell_activation_day0_rep2",
  "BCLLATLAS_14/Tcell_activation_day1_rep2"
)
t_act_list <- list()

for (lib in libraries) {
  # Read the data
  lib_path <- str_c("data/", lib, "/filtered_feature_bc_matrix/")
  expression_matrix <- readMM(str_c(lib_path, "matrix.mtx.gz"))
  barcodes <- read_csv(str_c(lib_path, "barcodes.tsv.gz"), col_names = FALSE)
  colnames(barcodes) <- "barcode"
  features <- read_tsv(str_c(lib_path, "features.tsv.gz"), col_names = FALSE)
  colnames(features) <- c("ensembl", "symbol", "feature_type")
  rownames(expression_matrix) <- features$symbol
  colnames(expression_matrix) <- barcodes$barcode
  
  # Separate HTO and RNA matrices
  hto_ind <- which(str_detect(features$feature_type, "Antibody Capture"))
  rna_ind <- which(str_detect(features$feature_type, "Gene Expression"))
  t_act_hto <- expression_matrix[hto_ind, ]
  t_act_rna <- expression_matrix[rna_ind, ]
  
  # Setup Seurat object
  t_act <- CreateSeuratObject(counts = t_act_rna)
  
  # Normalize RNA data with log normalization
  t_act <- NormalizeData(t_act)
  
  # Find and scale variable features
  t_act <- FindVariableFeatures(t_act, selection.method = "vst")
  t_act <- ScaleData(t_act, features = VariableFeatures(t_act))
  
  # Add HTO as an independent assay
  t_act[["HTO"]] <- CreateAssayObject(counts = t_act_hto)
  t_act <- NormalizeData(t_act, assay = "HTO", normalization.method = "CLR")
  
  # Demultiplex
  t_act <- HTODemux(t_act, assay = "HTO", positive.quantile = 0.99)
  
  # Append to list of Seurat objects
  t_act_list[[lib]] <- t_act
}
names(t_act_list) <- str_remove(names(t_act_list), "(BCLLATLAS_11/|BCLLATLAS_14/)")

5 Visualize

We can visualize the results as ridge plots or heatmaps:

ridge_l <- map(t_act_list, function(t_act) {
  Idents(t_act) <- "HTO_maxID"
  RidgePlot(
    t_act,
    assay = "HTO",
    features = rownames(t_act[["HTO"]])[1:6],
    ncol = 2
  )
})

heatmap_l <- map(t_act_list, function(t_act) {
  HTOHeatmap(t_act, assay = "HTO", ncells = 5000)
})

# Save
# walk2(ridge_l, names(ridge_l), function(ridge, lib) {
#   ggsave(
#     filename = str_c("results/plots/", lib, "_hashtag_demux_ridge.pdf"),
#     plot = ridge, height = 9, 
#     width = 16
#   )
# })
# walk2(heatmap_l, names(heatmap_l), function(heat, lib) {
#   pdf(file = str_c("results/plots/", lib, "_hashtag_demux_heatmap.pdf"), height = 4, width = 7)
#   print(heat)
#   dev.off()
# })
ridge_l
## $Tcell_activation_day0

## 
## $Tcell_activation_day2

## 
## $Tcell_activation_day0_rep2

## 
## $Tcell_activation_day1_rep2

heatmap_l
## $Tcell_activation_day0

## 
## $Tcell_activation_day2

## 
## $Tcell_activation_day0_rep2

## 
## $Tcell_activation_day1_rep2

6 Save Seurat objects

# saveRDS(t_act_list, "results/R_objects/t_act_Seurat_list_demultiplexed.rds")

7 Session Info

sessionInfo()
## R version 3.6.1 (2019-07-05)
## Platform: x86_64-apple-darwin15.6.0 (64-bit)
## Running under: macOS High Sierra 10.13.6
## 
## Matrix products: default
## BLAS:   /Library/Frameworks/R.framework/Versions/3.6/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/3.6/Resources/lib/libRlapack.dylib
## 
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
## 
## attached base packages:
## [1] parallel  stats4    stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
##  [1] forcats_0.5.0               stringr_1.4.0               dplyr_0.8.4                 purrr_0.3.3                 readr_1.3.1                 tidyr_1.0.2                 tibble_2.1.3                tidyverse_1.3.0             Matrix_1.2-18               Seurat_3.1.4                scater_1.14.6               ggplot2_3.3.0               SingleCellExperiment_1.8.0  SummarizedExperiment_1.16.1 DelayedArray_0.12.2         BiocParallel_1.20.1         matrixStats_0.55.0          Biobase_2.46.0              GenomicRanges_1.38.0        GenomeInfoDb_1.22.0         IRanges_2.20.2              S4Vectors_0.24.3            BiocGenerics_0.32.0         BiocStyle_2.14.4           
## 
## loaded via a namespace (and not attached):
##   [1] readxl_1.3.1             backports_1.1.5          sn_1.5-5                 plyr_1.8.6               igraph_1.2.4.2           lazyeval_0.2.2           splines_3.6.1            listenv_0.8.0            TH.data_1.0-10           digest_0.6.25            htmltools_0.4.0          magick_2.3               viridis_0.5.1            fansi_0.4.1              gdata_2.18.0             magrittr_1.5             cluster_2.1.0            ROCR_1.0-7               globals_0.12.5           modelr_0.1.6             RcppParallel_4.4.4       sandwich_2.5-1           colorspace_1.4-1         rvest_0.3.5              ggrepel_0.8.1            haven_2.2.0              xfun_0.12                crayon_1.3.4             RCurl_1.98-1.1           jsonlite_1.6.1           survival_3.1-8           zoo_1.8-7                ape_5.3                  glue_1.3.1               gtable_0.3.0             zlibbioc_1.32.0          XVector_0.26.0           leiden_0.3.3             BiocSingular_1.2.2       future.apply_1.4.0       scales_1.1.0             mvtnorm_1.1-0            DBI_1.1.0                bibtex_0.4.2.2           Rcpp_1.0.3               metap_1.3                plotrix_3.7-7           
##  [48] viridisLite_0.3.0        reticulate_1.14          rsvd_1.0.3               tsne_0.1-3               htmlwidgets_1.5.1        httr_1.4.1               gplots_3.0.3             RColorBrewer_1.1-2       TFisher_0.2.0            ica_1.0-2                farver_2.0.3             pkgconfig_2.0.3          uwot_0.1.5               dbplyr_1.4.2             labeling_0.3             tidyselect_1.0.0         rlang_0.4.5              reshape2_1.4.3           cellranger_1.1.0         munsell_0.5.0            tools_3.6.1              cli_2.0.2                generics_0.0.2           broom_0.5.5              ggridges_0.5.2           evaluate_0.14            yaml_2.2.1               npsurv_0.4-0             fs_1.3.2                 knitr_1.28               fitdistrplus_1.0-14      caTools_1.18.0           RANN_2.6.1               pbapply_1.4-2            future_1.16.0            nlme_3.1-145             xml2_1.2.2               rstudioapi_0.11          compiler_3.6.1           beeswarm_0.2.3           plotly_4.9.2             png_0.1-7                lsei_1.2-0               reprex_0.3.0             stringi_1.4.6            lattice_0.20-40          multtest_2.42.0         
##  [95] vctrs_0.2.3              mutoss_0.1-12            pillar_1.4.3             lifecycle_0.1.0          BiocManager_1.30.10      Rdpack_0.11-1            lmtest_0.9-37            RcppAnnoy_0.0.15         BiocNeighbors_1.4.2      data.table_1.12.8        cowplot_1.0.0            bitops_1.0-6             irlba_2.3.3              gbRd_0.4-11              patchwork_1.0.0          R6_2.4.1                 bookdown_0.18            KernSmooth_2.23-16       gridExtra_2.3            vipor_0.4.5              codetools_0.2-16         MASS_7.3-51.5            gtools_3.8.1             assertthat_0.2.1         withr_2.1.2              sctransform_0.2.1        mnormt_1.5-6             multcomp_1.4-12          GenomeInfoDbData_1.2.2   hms_0.5.3                grid_3.6.1               rmarkdown_2.1            DelayedMatrixStats_1.8.0 Rtsne_0.15               lubridate_1.7.4          numDeriv_2016.8-1.1      ggbeeswarm_0.6.0